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Abstract 

Synchronization between neuronal populations plays an important role in information transmission between brain areas. In 
particular, collective oscillations emerging from the synchronized activity of thousands of neurons can increase the 
functional connectivity between neural assemblies by coherently coordinating their phases. This synchrony of neuronal 
activity can take place within a cortical patch or between different cortical regions. While short-range interactions between 
neurons involve just a few milliseconds, communication through long-range projections between different regions could 
take up to tens of milliseconds. How these heterogeneous transmission delays affect communication between neuronal 
populations is not well known. To address this question, we have studied the dynamics of two bidirectionally delayed- 
coupled neuronal populations using conductance-based spiking models, examining how different synaptic delays give rise 
to in-phase/anti-phase transitions at particular frequencies within the gamma range, and how this behavior is related to the 
phase coherence between the two populations at different frequencies. We have used spectral analysis and information 
theory to quantify the information exchanged between the two networks. For different transmission delays between the 
two coupled populations, we analyze how the local field potential and multi-unit activity calculated from one population 
convey information in response to a set of external inputs applied to the other population. The results confirm that zero-lag 
synchronization maximizes information transmission, although out-of-phase synchronization allows for efficient 
communication provided the coupling delay, the phase lag between the populations, and the frequency of the oscillations 
are properly matched. 
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Introduction 

Brain function emerges from the collective dynamics of coupled 
neurons, the structural connectivity among which enables 
correlations between their firing activities. As a result of these 
correlations, effective neuronal networks function collectively on a 
mesoscopic scale, comprising thousands of coupled neurons that 
operate together, giving rise to emergent behavior. In awake 
animals, this collective dynamics takes the form of recurrent series 
of high and low network activity, corresponding with repetitive 
epochs of increased excitation over inhibition followed by boosts of 
inhibition. This leads to the appearance of rhythmicity at certain 
frequency bands. In particular, oscillations in the gamma-band 
(30Hz — 90Hz) are observed in several cortical areas in relation 
with cognitive tasks [1]. 

Synchronized oscillations can increase the functional connec- 
tivity between neural assemblies by coherently coordinating their 
firing dynamics. This hypothesis, known as communication 
through coherence (CTC), was proposed [2] as a mechanism by 
which gamma-band synchronization could regulate routing of 



information between brain areas. Since neuronal oscillations are 
associated with the dynamics of the excitatory-inhibitory balance, 
they represent periodic modulations of the excitability of neurons, 
which are more likely to spike within specific time windows (i.e. 
when inhibition is low). If two neuronal populations oscillate with 
a constant phase difference, then an effective transmission of 
information between the two groups of neurons is achieved 
provided the spikes sent by a population reach systematically the 
other population at the peaks of excitability. In that way, 
modulation of the relative phases of the emerging rhythms might 
dynamically generate functional cell assemblies [3-5] . 

A key requirement of the CTC mechanism is the existence of a 
constant phase difference between the two neuronal networks that 
reliably allows their binding, favoring communication. This 
coordination can be expected to arise from the synaptic coupling 
between the neurons of the two populations. But this coupling is 
not instantaneous, since propagation times between different 
cortical regions can reach up to several tens of milliseconds [6] . 
Previous CTC studies have mainly concentrated on describing the 
dependence of the coherence on the phase lag between the 
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Author Summary 

The correct operation of the brain requires a carefully 
orchestrated activity, which includes the establishment of 
synchronized behavior among multiple neuronal popula- 
tions. Synchronization of collective neuronal oscillations, in 
particular, has been suggested to mediate communication 
between brain areas, with the global oscillations acting as 
"information carriers" on which signals encoding specific 
stimuli or brain states are superimposed. But neuronal 
signals travel at finite speeds across the brain, thus leading 
to a wide range of delays in the coupling between 
neuronal populations. How the brain reaches the required 
level of coordination in the presence of such delays is still 
unclear. Here we approach this question in the case of two 
delay-coupled neuronal populations exhibiting collective 
oscillations in the gamma range. Our results show that 
effective communication can be reached even in the 
presence of relatively large delays between the popula- 
tions, which self-organize in either in-phase or anti-phase 
synchronized states. In those states the transmission 
delays, phase difference, and oscillation frequency match 
to allow for communication at a wide range of coupling 
delays between brain areas. 



neuronal populations [2,3,7], without examining systematically the 
relationship between the actual coupling delay and the phase lag at 
which the coherence is maximal. In fact, coupled nonlinear 
oscillators are known to adjust their phases upon frequency 
locking, leading under certain conditions to either in-phase (zero 
phase lag) or anti-phase synchronization (re-phase lag) [8]. Anti- 
phase patterns in cortical networks, for instance, have been widely 
studied [9]. Zero-lag synchronization, in turn, has been experi- 
mentally observed between gamma oscillations emerging from 
separated brain areas [10-12]. The conditions leading to zero-lag 
synchronization in neuronal oscillations are however somewhat 
stringent, requiring non- trivial spiking dynamics [13] or complex 
network architectures [14,15]. In particular, zero-lag synchroni- 
zation between two cortical areas has been shown to be possible 
even with long axonal delays [15,16], provided the two areas 
interact through a third oscillator, which could correspond to the 
thalamus [17,18]. 

But in contrast with most nonlinear oscillators neuronal 
populations are highly complex, since they contain a very large 
number of degrees of freedom (corresponding to the individual 
neurons), their oscillations are a pure collective phenomenon (the 
individual neurons in the population do not oscillate), and they 
operate in a broad frequency range. Additionally, neuronal 
populations are connected by a large number of axons, and 
inhomogeneities in the properties of those axons affect differen- 
tially the propagation speed of action potentials and lead to a wide 
spectrum of axonal delays rather than a uniform distribution [19]. 
It thus becomes necessary to study systematically the conditions 
under which two such complex oscillators synchronize (i.e. lock 
their frequencies), what is the resulting phase difference between 
them, how does this phase difference relate with the coupling delay 
(and with the frequency band being considered), and how is the 
efficiency of the communication between the two cortical areas 
affected by the delayed coupling. We address these questions in 
what follows. 

As mentioned above, within the CTC scenario effective 
communication arises when spikes from the emitting neuronal 
population reach the receiver population during the windows of 
maximum excitability. For this to happen two conditions have to be 



met: (1) the two coupled oscillators should be frequency locked, and 
(2) the transmission delay, the oscillation frequency, and the phase 
difference between the two oscillations should match. In particular, 
if the networks and the inter-connectivity is symmetric the second 
condition should hold in the two directions of spike propagation. 
The time delay (or rather, the distribution of time delays) is frxed as 
given by the anatomical connectivity. Therefore, it is the frequency 
of the oscillation spectrum what determines the particular phase lag 
that meets the matching condition. We have investigated whether 
this condition only occurs at specific rhythms, or if it holds at all 
frequencies. To this aim, we have represented mathematically two 
reciprocally connected identical neuronal populations using con- 
ductance-based models for both excitatory and inhibitory cells, and 
have studied how the heterogeneous axonal delays between the 
populations affect their synchronization. 

We have characterized the collective dynamics through a 
variable comparable to the local field potential (LFP) recordings 
[20]. In agreement with experimental data, the power of the 
modeled LFP decays with increasing frequencies [21]. Here we 
have focused on the particular dynamical regime in which the 
collective oscillations show a prominent contribution in the 
gamma range arising from the inhibitory (GABAergic) synaptic 
decay time constants [22]. Lower frequency bands contain a 
strong component arising from the noisy Poissonian distribution of 
interspike intervals (ISI), which affect the synaptic activation and 
hence do not reflect the intrinsic dynamics of the network. On the 
contrary, higher frequency bands of small amplitude reflect the 
fast dynamics of the action potentials, also affecting the synapse 
activation time course. 

The modeled neuronal networks exhibit other well-known 
features of cortical dynamics, such as coexistence of irregular firing 
at the single-neuron level with collective rhythmicity at the 
population level, arising from the synaptic recurrent connections 
between the excitatory and inhibitory neurons [23] (see Figure 1). 
The excitatory and inhibitory synaptic currents are balanced by 
compensating the higher number of excitatory neurons (80% of 
the whole network) with fast spiking inhibitory neurons and with 
strong inhibitory synaptic conductances. As a consequence, the 
neurons remain excitable but spent most of their time with a 
membrane voltage that fluctuates under the firing threshold. The 
gamma rhythm emerges from the periodic changes of this 
balanced synaptic current, which leads to periodic modulation of 
the distance to threshold. We have characterized the global 
activity of the network by means of averaging measures such as the 
aforementioned local field potential (LFP) and the multi-unit 
activity (MUA). We first used these measures to quantify phase 
coherence between the oscillatory activity of the two delay-coupled 
populations at varying mean axonal delays, observing transitions 
between in-phase and anti-phase dynamics. We next used 
information theory to quantify the response of one population 
(the receiver) to a varying external input impinging originally on 
the other population (the emitter). Our results show that 
information transmission is enhanced at zero-lag (in-phase) 
synchronization, and decreases at long delays for which commu- 
nication occurs through anti-phase dynamics. 

Results 

In-phase synchronization of collective oscillations under 
instantaneous coupling 

We start by considering an isolated population of 2000 neurons, 
of which 80% are excitatory and 20% are inhibitory. Each neuron 
forms on average 200 random connections within the network, 
and all pairs of coupled neurons exhibit a certain time delay, taken 
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Figure 1. Collective oscillations of a population of 2000 neurons. (A) Raster plot of 2000 neurons (in red the excitatory and in black the 
inhibitory neurons) for a 1500 — ms interval. (B) LFP time trace in a 1500— ms interval for an external mean rate of 7300 spikes/s. (C) MUA signal 
calculated counting the number of spikes of the neural population per unit time. (D) LFP power spectrum calculated using the Welch method 
averaged over 200 trials. The gray horizontal bar delimits the gamma peak band (30Hz — 52Hz). 
doi:10.1371/journal.pcbi.1003723.g001 



from a gamma distribution whose scale and shape parameters are 
both equal to unity. All neurons receive an external Poisson- 
distributed spike train whose instantaneous firing rate follows an 
Ornstein-Uhlenbeck process with a mean value set to 
7300spikes/s. This input and the excitatory recurrent synaptic 
activity are balanced by the recurrent inhibitory synaptic flow, 
since the GABAergic conductances are stronger than the 
glutamatergic AMPA ones. Furthermore, the inhibitory neurons 
fire at higher rates than the excitatory cells. Therefore, the 
membrane voltage of the neurons fluctuates below threshold, 
occasionally crossing it [24] . Despite the fact that the neurons fire 
sparsely and irregularly (see Figure 1A), a rhythmicity emerges 
when considering the dynamics of multiple action potentials 
elicited by thousands of neurons [23]. These oscillations represent 
the transient synchronized activity of neuronal assemblies, and can 
be revealed by population measures such as the local field 
potential (Figure IB) and the multi-unit activity (Figure 1C), 
defined in the Materials and Methods section. In the computa- 
tional model used throughout this work, the collective oscillatory 
dynamics exhibit a prominent gamma rhythm (Figure ID), whose 
period is mainly determined by the decay time constant of 
inhibition [23,25,26]. 

Another way of understanding the emergent gamma oscillations 
is by looking at the coupling between the MUA and the LFP. 
Since the LFP mainly captures the synaptic currents impinging on 
the pyramidal neurons (see Materials and Methods section), it is a 
measure of the excitability of the network. Hence, at those 
intervals in which inhibition is low (i.e. the inhibitory synaptic 
current fades away), the probability of firing is high. Due to the 



recurrent connections between the excitatory and inhibitory 
neurons, both the initiation and termination of the population 
bursts occur with a certain periodicity. Here this oscillatory pattern 
is around ~45Hz due to the inhibitory decay time constants [22]. 
The LFP and MUA are mutually locked to this frequency 
(Figure 2A), and the spikes occur with higher probability close to 
the troughs of the LFP (i.e. the minimum of inhibition, Figure 2B). 

We next consider two bidirectionally coupled neuronal 
networks of the type described above. Connections between the 
two areas are excitatory: 60% of the excitatory neurons of each 
network project randomly to 10% of the neurons belonging to the 
other pool. Although these parameter values cannot be general- 
ized to any two separate brain areas, for which the specific 
connectivity might determine their interactions, it is known that 
the probability of connection decays with distance [27-29]. Here 
we assume that the connectivity within a network is 2-fold the 
connectivity across networks, neglecting heterogeneity across 
neurons. Moreover, in order to obtain a certain amount of phase 
coherence between the two networks, we consider that the 
majority of excitatory neurons project onto the other network. A 
stronger (weaker) coupling will lead to unrealistically higher (lower) 
phase coherence values [30]. We have introduced time delays in 
the coupling between networks, assuming that the inter-areal 
delays are larger than the intra-areal delays due to long-range 
connections. We also consider that the inter-areal delays are 
distributed heterogenously across the system [19], following a 
gamma distribution whose mean and variance increase systemat- 
ically with the mean delay [15]. This accounts for the variability of 
transmission delays through axons with heterogeneous properties 
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Figure 2. Phase locking between LFP and MUA of a network. (A) LFP-MUA phase coherence for a single population. (B) Angle histogram of 

the phase difference between the LFP and MUA. The measures are averaged over 200 trials. 

doi:10.1371/journal.pcbi.1003723.g002 



(see Materials and Methods for the definition of the gamma 
distribution parameters). The mean inter-areal delay shown in the 
figures, hereafter termed T flY0 , accounts for the latency between the 
generation of a spike in a presynaptic neuron from one network 
and the elicitation of a postsynaptic potential in the other network. 

When coupled, the LFP power spectra of the two networks show 
the same gamma profile as in the absence of coupling, while the 
corresponding time series exhibit a substantial degree of correla- 
tion (Figure 3A inset). We next asked whether the broad spectrum 
of these neuronal oscillations allows for partial phase coherence to 
arise in specific frequency regions. Our phase coherence measure, 
described in the Materials and Methods section, quantifies 
between 0 and 1 the reliability of the phase difference A</> between 
pairs of oscillations, at a given frequency. Figure 3B shows the 
phase coherence between the LFPs of the two populations for 
instantaneous coupling (z axo = Oms). According to the regions of 
statistical significance observed experimentally [30], we considered 
phase coherence values above 0.08, which mainly occurs within 
the gamma band around the peak of the LFP power spectrum 
(horizontal gray bar in Figure 3B). This threshold corresponds to 
around four times the average phase coherence of the uncoupled 
case (see black dashed line in Figure 3B). 

We have also computed the time lag %j ag between the two 
signals (i.e. the time shift separating two equal phases of the 
coupled LFPs arising from each population) for all frequencies 
(Figure 3C), still in the case x axo = Oms. This time lag is only 
meaningful for significant phase coherence values that lead to a 
consistent A<j> across trials (red crosses in Figure 3B). The figure 
shows that for frequencies at which the phase coherence is 
significant, the LFP gamma rhythms of the two populations 
oscillate in phase (xi ag x0), i.e. the two LFPs are synchronized at 
zero lag. The error bars in Figure 3B,C represent the standard 
deviation across trials of phase coherence and x\ ag respectively, 



and are only shown for the region of significant phase coherence, 
since outside that region the phase distribution is very broad due to 
the variability across trials. Even within the significant region the 
standard deviation of Xj ag can be seen to decrease with increasing 
values of phase coherence, which confirms the inverse relation 
between phase coherence and the broadness of the phase 
distribution. 

Phase-coherence transitions for increasing coupling 
delay 

The fact that the two populations synchronize at zero lag when 
the coupling delay is zero is to be expected, and we now ask what 
happens in the presence of time delays. Figure 4 shows the phase 
coherence spectrum between the LFP oscillations for three 
different values of x axo . While phase coherence is again significant 
only around the gamma band (Figures 4A,C,E), the time traces 
look very different for small and large delays, with mostly in-phase 
dynamics for small delays (Figure 4B), whereas the populations are 
mostly in anti-phase for large delays (Figure 4F). For intermediate 
delays, interestingly, two coherence peaks appear (Figure 4C), and 
the corresponding time series exhibit both in-phase and anti-phase 
episodes (Figure 4D). These results indicate that in-phase 
dynamics seems to persist for non-zero coupling delays, eventually 
transitioning to an anti-phase regime with smaller, although still 
significant, phase coherence. Both types of dynamics seem to 
coexist for intermediate delays. 

In order to verify these conclusions, we have extended the 
analysis to a range of axonal delays, from Oms to 30ms, calculating 
the phase shift for the frequencies corresponding to both the peak 
of the power and the phase coherence spectra, termed f max - 
Figure 5A shows the value of the frequency at which the power 
spectrum is maximum, F,,, as a function of the coupling delay x axo . 
Note that varying x axo does not change the frequency peak of the 
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Figure 3. Collective oscillations of two coupled bidirectionally neural populations. The inter-areal axonal delay t„ A „ between the two 
neuronal pools is zero. (A) LFP time trace of the two populations in a 1000 — ms interval, for an external mean rate of 7300spikes/s. The inset shows 
the averaged time correlation of 200 LFP pairs. (B) Phase coherence between the LFPs of the two networks for varying frequency. The measure is 
averaged over 200 trials. The black dashed line represents the threshold (0.08) above which the phase coherence is considered significant (in red). (C) 
Time shift between the LFP oscillations of the networks for varying frequency. Red crosses show the time shifts corresponding to the frequencies at 
which the phase coherence is above threshold. The time shift is calculated as x/ ag = A(j>/27if„ mx , where A(j> is the phase difference at the frequency f max 
of maximum phase coherence. The gray bar delimits the gamma peak band (30Hz — 52Hz). The measure is averaged over 200 trials. 
doi:1 0.1 371 /journal.pcbi.1 003723.g003 



LFP power spectrum, which remains around 45 Hz for all 
coupling delays. We have added a gray bar delimiting the 
maximum power spectrum range within the gamma band 
corresponding to the extent of this local peak, highlighting the 
fact that the LFP gamma rhythm expands over a range of 
frequencies between approximately 30 — 52Hz. 

On the other hand, x axo clearly affects the frequency f max at 
which phase coherence is maximal, as shown by Figure 5B. In 
particular, f max exhibits a discontinuous jump around a coupling 
delay ~9ms, where two peaks of phase coherence coexist 
(consistent with the result shown in Figure 4C). The phase 
coherence values themselves are shown in color code in Figure 5C 
for different frequencies (vertical axis) and for varying z axo 
(horizontal axis). We have superimposed in that plot the line 
shown in panel A, which marks the maximum of the LFP power 
spectrum (black dashed line) within the gamma range, F y , as well 
as the whole extent of the peak (vertical gray bar). The local peaks 
of phase coherence fmax (black lines) corresponding to panel B are 
also superimposed to Figure 5C. 

For x axo = 0ms (as in Figure 3) the peak of phase coherence 
almost coincides with the peak of power spectrum. For increasing 
taxo, below 9ms, only the coherence peak at the lower frequency is 
significant, whereas between 10ms and 22ms only the coherence 
peak at the faster frequency is above threshold. The transition 



between these two regimes involves a coexistence of the local 
coherence peaks. We also observe that in both branches the 
frequency at which phase coherence is maximum f max decreases 
with the axonal delay, becoming clearly smaller than the gamma 
frequency peak Fy (dashed black line in Figure 5C). Making z axo 
greater than 22ms, which approximately matches the period of the 
power spectrum peak T y (l/F y «22ms), a new branch of phase 
coherence appears, thus leading again to coexistence of the two 
regimes. This emerging pattern is shown in Figure 5C for large 
inter-areal axonal delays and it is not marked in Figure 5B because 
the phase coherence is under the threshold. Hence, as z axo exceeds 
T y , the scenario of relative phases is repeated but now with cycle 
skipping. 

The phase coherence patterns shown in Figure 5C are affected 
by the inter-areal delay variability. If z axo is fixed to a constant 
value, the region of coexistence between the in-phase and anti- 
phase coherence patterns increases, and for delays approaching 
the oscillation period Ty the new peak emerging at Fy (detectable 
in Figure 5C and corresponding to in-phase dynamics in 
Figure 5E) becomes significant. This is shown in Supplementary 
Figure SIC, which displays the phase coherence for constant 
T axo = 20ms (blue line), in comparison with the case % axo = 0ms 
(black line) and the one with i axo drawn from a gamma 
distribution with mean 20ms (red line). 
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Figure 4. Phase coherence of two coupled bidirectionally neural populations for three different values of the inter-areal axonal 
delays x axo . Phase coherence spectrum and corresponding representative time series for T„. vo = 3ms (A,B), 9ms (C,D), and 17ms (E,F). The inter-areal 
delays follow a gamma distribution with a mean equal to corresponding inter-areal axonal delay r axo . The gray bars on the x-axes of plots A, C, and E 
delimit the gamma peak band (30Hz — 52Hz). The phase coherence measure is averaged over 200 trials. 
doi:10.1371/journal.pcbi.1003723.g004 



According to Figure 5C, maximum values of phase coherence 
fmax appear at different frequencies for each x axo . Significant 
values of phase coherence at a certain frequency can occur 
provided that there is a certain amount of spikes being 
simultaneously and reliably sent between the two networks. Since, 
by construction, the two neuronal pools are identical, the 
information flow can only be symmetrically transmitted for an 
in-phase, Aij) = 0, and/ or an anti-phase, A(j> = n, relationship 
between the two LFPs. Therefore, for any x axo we can obtain 
the corresponding f max that satisfies 2m axa f max = 0 or n. From this 
expression we can thus expect that larger x axo leads to smaller f max 
and that the anti-phase configuration is given by x axo equal to half 
the period corresponding to f max , not to be mistaken with T y /2, 
half the gamma period and equal to 1 1 ms. 

To verify the aforementioned remarks we have next calculated 

A(j) 

the time shift x/ ag between the two coupled LFPs as — — : . Figure 5D 

2nJ 

shows that, at the peak of the LFP power spectrum (here / = F Y ), 
Ziag is zero for low (0ms < x axo < 5ms) and large delays 
(1 7ms <x axo < 26ms). On the other hand, for intermediate 
(6ms < x axo < 16ms) and large delays (x ax0 >27ms) T/ a? corre- 
sponds to half the period of the gamma rhythm (T r /2=l/(2F r ) 



« 11ms). As mentioned before (see Figure 2), at frequency F r the 
MUA and the LFP in each population are frequency locked. 
Therefore, for any axonal delay, the presynaptic spikes arrive 
within the troughs of the postsynaptic LFP. We can interpret these 
sharp transitions from in-phase to anti-phase oscillations, appear- 
ing with a periodicity given by 7L as the way by which the system 
keeps the symmetry for any x axo . 

Since the maximum of phase coherence f max does not coincide 
with Fy, we have also obtained Xj ag along the peaks of phase 
coherence. Figure 5E confirms that only two patterns arise: in- 
phase and anti-phase, which can simultaneously occur in the 
region between 9ms and 10ms. The lowest frequency branch 
corresponds to T; B „ « 0ms, and thus to zero-lag synchronization. 
On the other hand, the highest frequency branch corresponds to a 
xi a g value that matches half the period of the corresponding 
frequency, i.e. 1 /(2f max ) (labeled by a red line in the plot), and thus 
corresponds to anti-phase synchronization. 

The full values of the time shift for all frequencies are shown in 
color code in Figure 5F. The region of zero-lag synchronization 
disappears as the delay increases, giving way to a region of anti- 
phase synchronization. Due to the oscillatory dynamics, for x axo 
greater than T y , frequencies close to the gamma peak are again 
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Figure 5. Phase coherence and time shift behavior in the case of bidirectional symmetric coupling for increasing inter-areal axonal 
delays x„ xo . (A) Frequency /"..(black arrow) at which the power spectrum is maximum and extent of the gamma peak (gray bar) (results for only one 
population are shown, since they are the same for both populations). (B) Frequencies at which the phase coherence exhibits local maxima, f max . (C) 
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compatible with an in-phase pattern. However, it is important to 
note that phase coherence is strongly decreased as the cycle is 
repeated again (x axo > T y ), probably due to loss of temporal self- 
coherence of the oscillations themselves. 

Thus, provided that the LFP-LFP phase coherence is significant, 
an effective coupling exists at which the two populations oscillate 
with a constant phase difference, which depends on the frequency 
and on the axonal delay. In particular, only two possible phase 
shifts are allowed, namely zero-lag (T; a „ x 0ms) and an anti-phase 
(xiag « 1 /(2/)) synchronization. 

Figure 5C shows that the frequency at which maximum phase 
coherence occurs, fmax, might not correspond to the predominant 
gamma rhythm at i 7 ] ,~45Hz, although it is close to it and within 
the extent of the gamma peak (gray vertical bar). Thus, phase 



coherence is bounded by the region in which spikes are still phase 
locked to the LFP (Figure 2). The separation between f max and F y 
is clear when x axo varies between 0 and T y /2. Phase coherence is 
achieved at slower rhythms that still reliably carry the action 
potentials. Hence, the spikes elicited by each population system- 
atically reach the other one at its excitability windows. Moreover, 
lower fmax implies larger excitability windows and allows the 
networks to be synchronized in phase. For larger x axo , corre- 
sponding slower frequencies lying outside the gamma peak do 
not efficiently transmit spikes, due to the bounded region in 
which MUA is locked to the LFP. Therefore, at large x axo 
the system moves towards an anti-phase configuration, where 
the time lag matches and compensates for the inter-areal axonal 
delay. 
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From phase coherence to communication 

The LFP oscillations studied so far are complex rhythms that 
convey a wide range of frequencies with a predominant 
component in the gamma range. We have seen before that the 
axonal delay x axo determines the relative dynamics of the coupled 
neuronal pools, which fall in either an in-phase or an anti-phase 
pattern. The phase relationship set by the two LFP signals is 
proposed to regulate the effectiveness of communication [2]. In 
other words, a stable phase difference determines the response 
of a neuronal population to inputs perturbing directly another 
area. Therefore, depending on the phase difference A(j> between 
two coherent LFPs, the response of the unperturbed population 
will replicate to a certain extent the response of the other 
population to the perturbation. We next study how, in the two 
different synchronization scenarios described above, inter-areal 
axonal delays affect information transmission during temporal 
windows, in which the phase difference and the frequency cannot 
be independent of each other. Note here the difference between 
phase coherence and effective communication. Rigorously speak- 
ing, communication occurs whenever spikes from one population 
arrive to the other one, and this is guaranteed provided that there 
is some coupling across networks. In contrast, effective commu- 
nication refers to a more specific situation in which information 
about the stimulus is being carried by the coupled populations. 

We can obtain a good proxy for communication by quantifying 
the response of a neuronal population (the receiver) to a 
perturbation that affects indirectly its dynamics via a second 
population coupled to it (the emitter), and which receives direcdy 
the perturbation. The success in communication can be observed 
in the transient amplification of the neuronal oscillations of the 
receiver [31]. The perturbation simulates different stimulus 
features, and consists of increases in the mean firing rate of the 
background synaptic activity impinging on a subpopulation of the 
emitter. We have examined, at different inter-areal axonal delays 
Ta.to, how well the LFP and MUA power spectra of the receiver 
convey information about the external stimulus being applied to 
the emitter. 

Since the connectivity within and between the two neuronal 
networks exhibits a certain degree of clustering, exciting a 
subpopulation of adjacent excitatory neurons from an area in 
the emitter population triggers a response in a well-defined 
subpopulation of neighboring neurons in the receiving population. 
We have chosen a set of different input intensities, 
5 = 8300,8800,9300,9800,10300,10800,11300 spikes/s, affecting 
400 long-range excitatory neurons from the emitter population 
during a 2-second time window. As a consequence of this 
perturbation, the amplitude of the LFP power spectrum increases 
with the strength of the perturbation (Supplementary Figures 
S2A,B, S3A,B and S4A,B), without altering the position of the 
gamma frequency peak (F y ~45Hz), consistent with the results 
were reported in [20]. 

Perturbing one of the populations breaks the symmetry of the 
system, since now the firing activity of the emitter is enhanced with 
respect to the receiver. As shown by the maps of phase coherence 
plotted in Figure 6, an increase of the external firing rate boosts 
phase coherence between the two LFPs. Moreover, the two 
frequency bands where phase coherence is significant merge into a 
single region at larger values of S concentrating closer to the 
gamma frequency peak F y ~ 45Hz. The corresponding xi ag values 
are shown in Figure 7 (note the different ranges of the axes, which 
now concentrate on the significant values of phase coherence to 
better observe the transition to the out-of-phase regime). 

At the gamma frequency peak F y the system undergoes a 
transition from in-phase to anti-phase dynamics as x axo increases. 



Small x axo lead to time shifts xi ag ~ 0 of the emitter's LFP relative 
to the receiver's LFP (Figure 7A-B) and thus, the two signals 
oscillate approximately in phase. However, the route to the anti- 
phase configuration changes as S is strengthened. In particular 
higher S trigger smoother transitions and the anti-phase regime 
becomes narrower. Figure 8B shows X/ ag values tracked at the 
gamma frequency peak F y = 45Hz. Here, larger S leads to a 
leader-laggard configuration in which the emitter LFP precedes 
the receiver LFP by a time lag that equals the axonal delay (see 
dashed black lines). Supplementary Figures S2C,D, S3C,D and 
S4C,D show the phase coherence and time shift for x axo = 3,9, and 
17ms (the same delays as Figure 4) for the whole range of 
frequencies. 

The dependence of the phase coherence on x axo for different S 
values is shown in Figure 6A-D, corresponding to a shift from a 
symmetric to an effectively asymmetric coupling. As the extra 
perturbation is applied only to one of the populations, the effective 
coupling approaches an unidirectional connectivity, although the 
structural links are not changed. This can be further explained by 
carrying on the same analysis in a structural unidirectional 
scenario, in which only one population projects onto neurons from 
the other network. Supplementary Figure S5A shows that 
increasing the delay i axo of the unidirectional transmission, the 
networks keep the phase difference constant at approximately the 
same frequency close to the power spectrum peak frequency. This 
represents a leader-laggard configuration and is similar to what 
happens in Figure 6D, where an over-excited subpopulation is 
driving the coupling between the two networks, still bidirectional 
but strongly asymmetric. The decrease of phase coherence with 
axonal delay is due to the variability in delay times: fixing x axo to a 
constant value of 20ms leads to maximal phase coherence values 
comparable to the case of no delay (Supplementary Figure SI A). 
Figure S5B shows that for increasing inter-areal axonal delays i axo , 
the time shift between the two synchronized networks (at 
frequencies corresponding to the significant phase coherence of 
Figure S5A) increases as long as x axo is smaller than half the period 
of LFP oscillation (l/f m a X ) and then approaches zero, thus leading 
again to a transition from in-phase to anti-phase synchronization 
at frequencies close to that of the power spectrum peak F y . 

Phase coherence can influence the transmission of information 
between neuronal populations. As mentioned in the Introduction, 
the CTC hypothesis [2] suggests that phase relations between 
coupled areas modulate the response of a receiver area to 
presynaptic input coming from an emitter area. In order to 
maximize this response, the axonal delay x axo , the frequency / of 
the oscillations and the phase difference A<j> should verify 
A(j> = 2nx axo f. When this relationship holds, spikes fired in the 
emitting population at a specific phase of the signal (for instance at 
the troughs of the LFP, which correspond to the maxima of 
excitability) arrive at the receiving area at the same phase (and 
thus at the same excitability maximum), triggering a maximal 
response in the receiving area. On the contrary, if A(f> does not 
fulfill the relationship given above (or if it randomly varies), 
effective communication will not be achieved [31]. This condition 
is relevant at the frequencies at which the MUA and the LFP are 
phase locked (Figure 2). Otherwise, the troughs of the LFP do not 
correspond to intervals of maximum firing within the same 
population, and the peaks of MUA do not occur reliably with the 
same periodicity as the LFP. 

In order to quantify the efficiency of communication, we have 
computed the mutual information I(S; Rf) (defined in the 
Material and Methods section) between the power spectrum Rf 
at a frequency / of both the LFP and MUA of the receiver and the 
set of stimuli S applied to the emitter. This definition of 
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Figure 6. Phase coherence in the case of bidirectional asymmetric coupling for increasing extra inputs. Phase coherence between LFPs 
of the two networks, in color code, as a function of frequency (y-axis) and of the inter-areal axonal delay r„ v „ (x-axis) for different stimuli: (A) 8300, (B) 
9300, (C) 10300, (D) 11300 spikes/s. The measures are averaged over 200 trials for each i m and stimulus. 
doi:1 0.1 371 /journal.pcbi.1 003723.g006 



information does not require any assumption about the stimulus 
features being encoded by the neural signals [32,33]. I(S;Rf) 
quantifies the reduction of the uncertainty in predicting the 
applied stimulus given a single observation of the triggered 
response, and is measured in units of bits (1 bit means a reduction 
of uncertainty of a factor of two). Several experiments have been 
done with this tool to characterize information transmission in the 
primary visual cortex of macaques in response to a naturalistic 
stimulus [33]. Several other studies have been performed using the 
LFP power spectrum as a measure of mutual information, showing 
the usefulness of this approach both experimentally and compu- 
tationally [20] . The advantages of this approach are described in 
detail in [34,35]. 

To compute I(S; Rf), we have run extensive simulations to 
properly estimate the conditional probabilities used in mutual 
information measures. The techniques adopted in order to reduce 
the bias error of the estimation of conditional probability due to 
the finite number of samples are explained in the Material and 
Methods section. Figure 9 shows that the mutual information is 
non-negligible only within the gamma range (p<0.05; bootstrap 
test), in a narrow region around the peak of the power spectrum 
F v . This is consistent with the fact that the emitter encodes the 



different stimulus strengths in the gamma band, i.e. other regions 
of the LFP power spectrum are not affected (Figure S2-S4A). 
Therefore, information transmission occurs within the gamma 
peak (the mutual information spectrum of the two networks, 
computed from the LFP, for z axo = 3,9,17ms is plotted in 
Supplementary Figures S2E-S4E). Moreover, functional interac- 
tions between coupled neuronal populations can be maximized 
if their phase difference is close to zero, i.e. for short axonal 
delays. 

While I(S; Rf) is lower when computed for the LFP power 
spectrum (Figure 9A) than for the MUA power spectrum 
(Figure 9B), it decreases monotonically in both cases for increasing 
axonal delays. This behavior contrasts with the one shown in 
Figure 5C, in which the maximum phase coherence in the absence 
of stimulus occurs at varying frequencies f max for different Zaxo- 
Moreover, f max lies outside the significant mutual information 
spectrum. However, at large enough S the phase coherence 
pattern (Figure 6D) closely resembles the mutual information 
dependency with z axo (Figure 9), since here f max = F y . 

Mutual information encoded in the power spectrum is bounded 
to the frequencies at which spikes are maximally frequency locked 
(Figure 2). Although this measure does not take into account the 
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Figure 7. Time shift in the case of bidirectional asymmetric coupling for increasing extra inputs. Effective time shift in milliseconds 
between LFPs of the two networks, in color code, as a function of frequency (y-axis) and of the inter-areal axonal delay t„ ao (x-axis) for different 
stimuli: (A) 8300, (B) 9300, (C) 10300, (D) 11300 spikes/s. The measures are averaged over 200 trials for each z m and stimulus. 
doi:10.1371/journal.pcbi.1003723.g007 



phase difference between the two LFP signals, their dynamics 
clearly rely on their relative time lag. Therefore, significant phase 
coherence is needed in order to reliably connect in time the 
excitability time windows of both networks, but it is not sufficient 
to achieve a maximal response of the receiver. In order to meet 
this second requirement, the frequency at which phase coherence 
is obtained needs to carry a precise timing of the action potentials, 
otherwise the presynaptic current will not elicit a postsynaptic 
response. Even the emitter population can only encode the 
stimulus strength as variations in the amplitude of the gamma 
frequency peak, since it is at Fy that changes in the LFP represent 
synchronized alterations in the MUA. 

A symmetric coupling scenario allows for two emerging stable 
regimes, in-phase A^/„ g = 0 and anti-phase A(j) = n, while in an 
asymmetric regime the most excitable network leads the dynamics 
(Tlag = T a xo)- Therefore, in the presence of axonal delays, the latter 
case is not compatible with the in-phase/anti-phase condition. 
The symmetry breaking allows for ii ag to follow z axo , increasing 
phase coherence at the gamma rhythm and thus the receiver's 
response. In summary, efficient communication needs a sufficient 
locking between the spikes being transmitted and the LFP, the 



carrier of this information. This is maximized at the gamma 
frequency peak Fy, here ~45Hz, at which changes in the power 
spectrum due to external stimuli become particularly evident. The 
coupling axonal delay % axo modulates the level of phase coherence 
within all the gamma range, and strong driving of one of the 
populations precisely favors the ~45Hz frequency channel. As 
observed above, the variability of axonal delay affects the regions 
where the phase coherence maximum is significant. Supplemen- 
tary Figures S6A,B show the LFP and MUA mutual information 
in the unidirectional case. As in the case of bidirectional coupling, 
the flow of information occurs at Fy, where the MUA and LFP are 
frequency locked and the emitter encodes the stimulus strength. 
Specially, mutual information is higher at small z axo , when the 
networks are synchronized in phase. In the unidirectional 
configuration the mutual information shows a strong dependence 
on Xaxo, as in the case of phase coherence discussed above. This is 
due again to the variability of axonal delays. For a fixed time 
delay, the mutual information in the unidirectional coupling case 
does not show a substantial decrease for increasing % axo 
(Supplementary Figure SIB). The bidirectional configuration also 
exhibits a less significant decrease of the mutual information 
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Figure 8. Time shift behavior at the peak of power spectrum for increasing inter-areal axonal delays for different extra inputs. Effect 
of the external input perturbation on the coupled neuronal populations for increasing stimulus strengths 5 = 8300,8800,9300,9800,10300,10800, 
1 1300 spikes/s. (A) Frequency of the power spectrum peak. (B) Time shift corresponding to spectral peak frequency. The dashed lines show the ideal 
cases for which t] a g=t axo and its anti-phase equivalent. 
doi:10.1371/journal.pcbi.1003723.g008 



maximum for constant increasing x axo (Supplementary Figure 
SID). This is consistent with the phase coherence peak 
corresponding to in-phase dynamics that persists for increasing 
constant delay (Supplementary Figure SIC). 

Our results show that phase coherence cannot be taken as a 
precursor of information transmission. Phase coherence can be 
achieved in a broad range of frequencies around the gamma peak 
F y (Figure 6). Therefore, the spikes impinging on each network are 
able to bound the two populations in a constant phase relationship, 



constrained by the symmetry of the effective coupling. However, in 
order to communicate, presynaptic spikes must trigger a postsyn- 
aptic response. This requires that the presynaptic action potentials 
are synchronized in time to facilitate the integration of the synaptic 
currents. Hence, changes in the LFP and MUA amplitude occur 
precisely at F y and mutual information does the same (Figure 9). 
Stimulus that are able to modify the response of a population within 
a wider frequency range (i.e. not frequency specific) could possibly 
alter the situation here described. 
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Figure 9. Mutual information carried by LFP and MUA power spectrum of the receiver. Mutual information between the set of stimuli 
5 = 8300,8800,9300,9800,10300,10800,11300 spikes/s and the neural response given by the LFP (A) and MUA (B) power spectra for increasing 
coupling delays x axo . The gray arrow in the color scale refers to significance threshold (/7<0.05, bootstrap test). The measures are averaged over 200 
trials for each x„ xo and stimulus. 
doi:1 0.1 371 /journal.pcbi.1 003723.g009 
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Discussion 

Here we have examined how heterogeneous inter-areal 
synaptic delay influences the coupling between the collective 
dynamics of two neuronal populations. To that end, we first used 
population measures such as the local field potential and the 
multi-unit activity, by analogy with experimental studies, to 
capture the collective oscillatory dynamics of individual neuronal 
populations. In the presence of excitatory coupling, the LFP and 
MUA activities of two identical delayed neuronal networks 
oscillate in the gamma range, with a significant broad peak 
between 30 and 52Hz, which does not depend on the axonal 
delay x axo - The emergence of this gamma peak in the isolated 
populations is due to the recurrence between excitatory and 
inhibitory synaptic activity, as revealed by the phase locking 
between the LFP and MUA signals (Figure 2). In contrast with 
the power spectrum, phase coherence is strongly affected by the 
axonal delays between the populations. We have seen that in- 
phase and anti-phase patterns occur at various frequencies for 
different ranges of x axo , with high values of phase coherence 
occurring at frequencies below the gamma frequency peak F y 
(Figure 5). 

The phase coherence pattern shown in Figure 5C corresponds 
to a pure symmetrical connectivity, in which both the structural 
and functional coupling are equal in both directions (in contrast 
with the unidirectional case of Figure S5). The reciprocity between 
the feedback and feedforward pathways across cortical areas is not 
an unrealistic assumption [36], although the specificity of such 
synapses might differ in each direction in order to account for the 
different effects of the top-down and bottom-up projections. Here 
we show that increasing axonal delays x axo lead to either an in- 
phase or anti-phase synchronization with a vanishing maximal 
phase coherence at frequencies f max below F y although lying 
within the gamma peak. Hence, in basal conditions, there is always 
a certain reliable phase relationship, provided x axo is small, relative 
to the period T y . 

The interesting point raised by the communication through 
coherence hypothesis [2] , is whether phase coherence can forecast 
efficient communication between two populations in the presence 
of a stimulus. According to the modulatory role of the top-down 
pathway, attention can determine which visual cues we are aware 
of [37,38]. In principle two situations can arise: either a stimulus 
catches our attention (such as an unexpected noise or object) or we 
are being attentive to an expected stimulus (such as waiting the 
traffic light to turn green). In the first situation the communication 
outline between a primary cortical area and the associative areas is 
driven by the stimulus, while in the second case it is due to the 
internal cognitive state. The firing activity in visual areas has been 
shown to significantly increase even in the attentive state prior to 
the stimulus presentation [39]. Hence our results, in which we 
progressively increase the firing rate impinging on one population, 
could be viewed as arising from the alteration of the underlying 
attentional state. 

The experimental study of [38] shows that a neuronal cell 
assembly in V4 is phase coherent with an area in VI that responds 
to a selected stimulus, while it is not with a VI area that is not 
driven by the input. Here we have not modeled a competitive 
scenario between two networks. Instead we have focused on the 
mechanisms by which two neural pools can modulate their 
communication when they are simultaneously oscillating in the 
gamma band. We have quantified the efficiency of communication 
between the two neuronal networks as the ability of a population 
to encode information of an input which perturbs directly another 
coupled population. Mutual information measures between either 



the LFP or MUA power spectrum and the set of applied stimuli S 
show that significant values concentrate around the gamma 
frequency peak (« 45Hz). Mutual information decreases for long 
inter-areal axonal delays, and is slightly lower when the neural 
response is characterized by the LFP power spectrum than by the 
MUA power spectrum. 

Despite the fact that the LFP reflects the afferent and local 
synaptic currents within a given neuronal network, and that the 
MUA only captures the action potentials within this network, these 
two signals are closely related. As mentioned above, the gamma 
LFP rhythm reflects the dynamics of the excitatory balance. 
Increases in inhibition silence the spiking activity and therefore the 
MUA signal, although the GABAergic current is enhanced. 
Decreases in inhibition boost the spiking activity and therefore the 
MUA signal, although the GABAergic current is reduced. The 
peak at 45Hz in the LFP-MUA phase coherence (Figure 2) reveals 
this phase locking between the two signals. 

The arrival of each set of presynaptic spikes perturbs the 
postsynaptic LFP and might or might not elicit a postsynaptic 
suprathreshold response (captured by the postsynaptic MUA) 
depending on the mean distance to the excitatory threshold of the 
populations. Bursts of activity occur at each pool with a 
periodicity that fluctuates within the gamma band and are locked 
to the troughs of the LFP at this frequency. According to the 
CTC hypothesis, maximum communication requires that spikes 
from each population reach the peaks of excitability of the target 
area simultaneously in both coupling directions. Thus, efficient 
communication is restricted to the gamma peak, as revealed by 
the mutual information (Figure 9) and preferentially at relatively 
small x axo . This condition is only met for in -phase or anti-phase 
synchronization of the gamma rhythm: small axonal delays x axo 
are able to tie two LFP troughs only at zero-lag synchronization, 
while larger T axo require anti-phase synchronization. In principle, 
as i axo increases zero-lag synchronization could again mediate 
communication by skipping one cycle. However, due to loss of 
phase consistence, mutual information decays with increasing 

Xa X o • 

Here we show that phase coherence emerges spontaneously due 
to the excitatory coupling between areas without the need of 
further constrains (Figure 5C). Higher stimulation of a particular 
population (the emitter), which enhances the LFP power spectrum 
amplitude of the gamma peak, increases the range of phase 
coherence to larger axonal delays (Figure 6D). The delay 
determines the phase shift between the two signals, with the 
emitter leading the oscillations. According to [38], phase 
coherence is revealing communication in the sense of spike 
propagation, which in our case extends to frequencies within the 
gamma peak. However, efficient communication in the sense of the 
information encoding in the postsynaptic response, is restricted to 
a narrower band (Figure 9) that maximizes spike synchronization. 
Adopting a spectrum of delays with increasing variability for 
increasing values of x axo , instead of an (unrealistic) constant delay, 
affects quantitatively the results of phase coherence and mutual 
information but does not produce any strong qualitative change in 
the findings of the paper. However the effect of variability cannot 
be ignored, given the dispersion of axonal delays observed in 
experimental studies [19]. 

Figure 10 shows a schematic diagram of the two oscillatory 
LFPs filtered around the gamma frequency peak 
(l/Fy = Ty = 22 ms) with the bursts of spikes locked at their 
troughs in agreement with Figure 2. For a delayed coupling, zero- 
lag synchronization does not lead to a symmetric configuration 
demanding that the two oscillations are reciprocally influenced at 
the same phase. Therefore the system rapidly shifts toward an 
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Figure 10. Carriers of information and signals. Diagram of two oscillatory LFPs filtered around the power spectrum peak («45Hz), with a short 
spike train locked at their troughs for different r am : (A) x„ xo = 3ms, representing zero-lag synchronization and (B) z a v „ = 9ms, representing anti-phase 
synchronization. 

doi:10.1371/journal.pcbi.1003723.g010 



anti-phase synchronization at which z ax0 roughly equals half of 
the period of the LFP (Figure 10B). Importantly, when the 
symmetry of the system is broken (for instance by perturbing just 
one of the populations), the possible stable solutions are no longer 
the in-phase or the anti-phase regime. In this case, phase 
coherence can be achieved through a leader-laggard configura- 
tion in which the time lag equals the inter-axonal delay. Without 
the symmetry constraint, this situation is achieved at the gamma 
frequency peak, for which the spikes of each population are 
preferentially locked to the LFP and changes in their power 
spectrum are maximized. 

In conclusion, we have studied two neuronal populations 
coupled synaptically with non-negligible delays. Our modeling 
results show that the populations organize their joint collective 
dynamics in patterns of in-phase or anti-phase synchronization, 
depending on the delay. Unidirectional couplings, either structural 
or functional, lead to a leader-laggard configuration with an out- 
of-phase synchronization determined by the axonal delay. Our 
study shows the dichotomy between phase coherence and 
communication. Whereas phase coherence arises due to LFP 
phase perturbations through the propagated spikes, communica- 
tion is caused here by an increase in the firing response. The first 
occurs at different frequencies for every z axo in order to conserve 
the functional connectivity. The second requires the spikes to be 
tightly locked to the LFP and at a faster frequency F y to enable 
spike integration, and hence a signal response that can be 
synaptically propagated. 



Materials and Methods 

Computational model 

We consider two populations of 2000 neurons, 80% of which 
are excitatory while the remaining 20% are inhibitory [40] . Each 
neuron connects on average with 200 other cells through only 
chemical synapses. The structural connectivity is built according 
with the Watts-Strogatz small-world algorithm [41]. The rewiring 
probability is set to 0.5, so that the connectivity shows a certain 
degree of clustering, which favors the connections between 
neighboring neurons. Coupling between the two networks is 
mediated by 60% of the neurons of one population making 
random long-range excitatory projections with 10% of the neurons 
belonging to the other population. Here we assume that the 
connectivity within a network is 2-fold the connectivity across 
networks, neglecting heterogeneity across neurons. Moreover, in 
order to obtain a certain amount of phase coherence between the 
two networks, we consider that the majority of excitatory neurons 
project onto the other network. A stronger (weaker) coupling will 
lead to unrealistically higher (lower) phase coherence values [30]. 
We introduced a synaptic transmission delay within and among 
the networks, taken from a gamma distribution, assuming that 
internal delays (taken from a gamma distribution whose scale and 
shape parameters are fixed to 1) are smaller than the inter-area 
delays. The axonal delays, termed x axo in the paper, stand for the 
time between the generation of a spike in a presynaptic neuron 
from one network and the elicitation of a postsynaptic potential in 
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the other network. These delays are taken from a gamma 
distribution whose mean and variance increase with increasing 
x axo- We choose the scale parameter of the distribution equal to 
unity, so that the shape parameter equals i axo . In this way the 
coefficient of variation (CV) decreases for increasing mean as 
l/y / niean. In our analysis we systematically vary x axo between 
Oms and 30ms. 

Each neuron is dynamically described by the Hodgkin and 
Huxley (HH) model. The dynamics of the membrane voltage is 
given by: 

dV , , 

C—^ = -g Na m i h(V-E Na )-g k if{V-E k ) + ,^ 

-g L (V-E L )+ I ext + I syn , 

where C = 0.25nF (0.50nF) is the membrane capacitance for 
inhibitory (excitatory) neurons, the constants gff a = 12.5/iS, 
gx = 4.74/iS, and gx = 0.025^S are the maximal conductances 
of the sodium, potassium, and leakage channels, respectively, and 
E Na =40mV, E K =-S0mV, and E L = -65mV stand for the 
corresponding reversal potentials. According to the HH formula- 
tion, the voltage-gated ion channels are described by the following 
set of differential equations 



dm 
dt 



dh 
dt 



-a m {V)(\-m)-p m (V)m, 



= z h (V)(\-h)-ft h {V)h, 



(2) 



Given that m activates rapidly, we replace it by its steady-state 

value Tflco = — 7T- 

0im+P m 

In Equation (1) I syn is the synaptic current coming from the 
neighboring neurons impinging on a neuronal cell: 



Isyn(t)=g S yn{t)(V{t)-E syn ) 



(4) 



where g syn {i) is the synaptic conductance and E syn is the reversal 
potential of the synapse. For positive values of E syn the synapse is 
depolarizing or excitatory (E syn = OmV for glutamate receptors), 
otherwise it is hyperpolarizing or inhibitory [E syn = — 70mV for 
GABA receptors). In the equation (4) the synaptic conductance is 
described by: 



gsynif) 



1 decay ' 



7 Z decay 



-<->j 
-C x rise 



(5) 



where Tdecay an d 1,-ise ar e the decay and rise synaptic time, 
respectively, and gsyn is tuned in order to obtain a balance between 
excitation and inhibition. The constant g syn is set to maintain the 
postsynaptic potential (PSP) amplitudes within physiological 
ranges. All parameters values can be found in [26,42] . 

In equation (1) I ext represents an heterogenous Poisson train of 
excitatory presynaptic potentials with a mean event rate that varies 
following an Ornstein-Uhlenbeck process (see Supplementary 
Figure S7). This incoming external current mimics the direct input 
from any other area external to the network considered here. The 
instantaneous rate, X{i), of the external excitatory train of spikes is 
generated according to an Ornstein-Uhlenbeck process, as 
considered in [20]: 



-=z h (V)(\-ri)-P n (V)n, 



where the gating variables m(t), h(i) and n(t) represent the 
activation and inactivation of the sodium channels and the 
activation of the potassium channels, respectively. The voltage- 
dependent transition rates are 



^ F) = °- 1 l-exp(-(K + 16)/10) ' (3) 



A„(F) = 4exp(-(F+41)/18), 



a h ( V) = 0.07 exp ( - ( V + 30) /20), 



^ = -A(0 + a(0yf?(0 (6) 

where a(t) is the standard deviation of the noisy process and is set 
to 0.6 spikes/s. x is set to 16ms, leading to a power spectrum for 
the X time series that is flat up to a cut-off frequency 

/= =9.9Hz. t](t) is a Gaussian white noise. 

ZTtT 

The model has been integrated using the Heun algorithm [43], 
with a time step of 0.05ms. All simulations represent 2.0 seconds of 
activity. The connectivity, initial conditions and noise realization 
were varied from trial to trial. 

LFP and MUA 

We quantified the activity of the network in two different ways. 
We calculated the multi-unit activity (MUA) as the total number of 
spikes per unit time in the population, and the local field potential 
(LFP) as the sum of the absolute values of the excitatory and 
inhibitory synaptic currents acting upon the excitatory neurons, 
averaged over this population [20]: 



ft,(F) = [l+exp(-F710)]- 



LFP=R e (\I AMPA \ + \I GABA \), 



(7) 



a„(K) = 0.01 



(F + 20) 



l-exp(-(F+20)/10)' 



P„( V) = 0. 125 exp ( - ( V + 30)/80). 



where <...> denotes the average over all excitatory neurons. The 
term Iampa accounts for both the external excitatory heteroge- 
neous Poisson spike train and the recurrent excitatory synaptic 
current due to the network, while Igaba corresponds to the 
recurrent inhibitory synaptic current. R e represents the resistance 
of a typical electrode used for extracellular measurements, here 
chosen to be 1 MQ. 
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Spectral analysis 

We computed the power spectral density of LFPs and MUAs 
using the Welch method: the signal is split up into 256 point 
segments with 50% overlap. The overlapping segments are 
windowed with a Hamming window. The modified periodogram 
is calculated by computing the discrete Fourier Transform, and 
then computing the square magnitude of the result. The modified 
periodograms are then averaged to obtain the PSD estimate, 
which reduces the variance of the individual power measurements. 
The code has been implemented in MATLAB. Spectral quantities 
are averaged over 200 trials and phase coherence over 1000 trials. 

Phase coherence 

Phase coherence is calculated as in [30]: 



C xy (f)- 



1 S xy (f,n) 
N^\S xy (f,n)\ 



(8) 



where x and y denote the two signals, and S x) ,(f,n) is the cross- 
spectrum between them. Since in each trial the cross spectral 
density is normalized by its amplitude, each term of the sum is a 
unit-length vector representation of the phase relation A(j>(f,n). In 
other words, A<f>(f ,n) = <l> y — <j> x is the phase lag between the two 
signals at frequency / in the data segment n. Hence C xy (f) 
quantifies how broad is the distribution of A<j>(f,n) within the 2n- 
cycle. Averaging A(j)(f,ri) across all N data segments provides a 
mean angle A<f>(f ). In our work A<j>(f ) is converted into a time 
shift, termed Ti ag in the paper, dividing by the corresponding 



frequency t(f)~ 



2nf ■ 



This quantity measures the time separa- 



tion between an LFP maximum in one population and the 
following maximum belonging to the other population. 

Mutual information 

An important mathematical tool to quantify information 
transmission in noisy systems is provided by information theory. 
We calculate the Mutual Information I(S; R) between the 
stimulus S and the response R as follows. The broadband LFP 
signal reproduces the variations in neural population activity over 
a wide range of time scales [44] . Thus LFPs signals are useful to 
qualitatively characterize mechanisms of information processing, 
because it is possible through them to verify if there are priviliged 
time scales for information processing. We can think that 
information is spread over all frequencies, or that each frequency 
contributes separately to the information representation. Given 
that we are interested in how the collective dynamics of the 
population carries information, we quantify the neural response 
Rf as the power of either the LFP or the MUA at frequency f, and 
we consider as stimuli different external firing rates impinging on 
one of the two populations. We compute the information between 
the stimulus S and the response Rf as: 



7(S; R f )=Y t P(^P(r f -|j)log : 



7 



P(r f \s) 
P{r f ) 



(9) 



where P(s) is the probability of presenting stimulus s (equal to the 
inverse of the total number of different external firing rates, 
namely of stimuli), P(jf) is the probability of observing power rf 
across all trials in response to any stimulus, and P(rf\s) is the 
probability of observing power rf at frequency / in response to a 
single stimulus S. I(S; Rf) quantifies the reduction of uncertainty 



about the stimulus that can be gained from observing a single-trial 
neural response, and we measured it in units of bits (1 bit means a 
reduction of uncertainty of a factor of two) [35] . This measure 
allows us to evaluate how well the power Rf of either the LFP or 
MUA encodes the stimulus at a certain frequency /. 

To facilitate the sampling of response probabilities, the space of 
power values at each frequency was binned into 6 equipopulated 
bins [33]. We used seven different firing rates of the external 
Poisson-distributed input, for a time T = 2s. An important issue to 
be solved regarding the calculation of the theoretical mutual 
information is that it requires knowledge of the full stimulus- 
response probability distributions, and obviously these probabilities 
are calculated from a finite number of stimulus-response trials. This 
leads to the so-called limited sampling bias, which corresponds to a 
systematic error in the estimate of information. We used the method 
described in [45] to estimate the bias of the information quantity 
and then we checked for the residual bias by applying a bootstrap 
procedure in which mutual information is calculated when the 
stimuli and responses are paired at random. If the information 
quantity is not zero (it should be in the case of non finite samples), 
this is an indication of the bias and the bootstrap estimate of this 
error should be removed from the mutual information. After 
applying these procedures, the information quantity estimation 
could be defined as significant. Several toolboxes provide different 
bias-correction techniques, which allow accurate estimates of 
information theoretic quantities from realistically collectable 
amounts of data [46,47]. In order to accomplish those tasks, we 
used the Information Breakdown Toolbox (ibTB), a MATLAB 
toolbox implementing several information estimates and bias 
corrections. It does this via a novel algorithm to minimize the 
number of operations required during the direct entropy estimation, 
which results in extremely high speed of computation. It contains a 
number of algorithms which have been thoroughly tested and 
exemplified not only on spike train data but also on data from 
analogue brain signals such as LFPs and EEGs [47]. 

Supporting Information 

Figure SI Phase coherence for constant inter-areal 
delay, (a) Phase coherence between the two LFP oscillations in 
the unidirectional coupling configuration when x axo = 0ms (black 
line), T axo = 20ms (blue line) and % axo is taken from a gamma 
distribution of mean 20ms (red line), (b) Mutual information 
between the set of stimuli 5 = 8300,8800,9300,9800,10300,10800, 
11300 spikes/s and the neural response given by the LFP in the 
same unidirectional coupling configuration, (c) Phase coherence 
and (d) mutual information in the bidirectional coupling 
configuration. Phase coherence measures are averaged over 
1000 trials. Mutual Information measures are averaged over 5 
sets of 200 trials for each stimulus. 
(EPS) 

Figure S2 Effect of external stimulation for small 
coupling delay, (a) LFP power spectrum of the directly 
stimulated population for different external rates (8300, 
8800,9300,9800,10300,10800,11300 spikes/s). (b) LFP power 
spectrum of the second population, (c) Phase coherence between 
the two LFPs for different external rates, (d) Effective time delay 
between the two pairs of LFP oscillations at frequencies where the 
phase coherence is significant for different external rates, (e) 
Mutual information between the LFPs of the two populations. Red 
dashed line corresponds to significance threshold (/?<0.05; 
bootstrap test) for information. The mean inter-delay between 
the pools is « 3ms. The measures are averaged over 200 trials. 
(EPS) 
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Figure S3 Effect of external stimulation for intermedi- 
ate coupling delay. The meaning of the plots is the same as in 
Suppl. Fig. S2. The mean inter-delay between the pools is here 
x9ms. 
(EPS) 

Figure S4 Effect of external stimulation for large 
coupling delay. The meaning of the plots is the same as in 
Suppl. Fig. S2. The mean inter-delay between the pools is here 
« 17ms. 

(EPS) 

Figure S5 Phase coherence and time shift in the case of 
unidirectional coupling, (a) Phase coherence, in color code, as 
a function of frequency (y-axis) and of the inter-areal axonal delay 
taxo (x-axis) in the case of unidirectional coupling from the emitter 
to the receiver, (b) Time shift Tlag, in color code, as a function of 
frequency (y-axis) and of the inter-areal axonal delay z axo (x-axis) 
in the case of unidirectional coupling from the emitter to the 
receiver. The measures are averaged over 1000 trials consistently 
with the symmetric coupling. 
(EPS) 
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